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and  optimization  of  core  shaping  to  provide  specific  eddy-current  flow 
patterns  Is  given. 


Accession  For 

MTIS  GRAtl 

/ 

dtic  tab 

□ 

Uivinnounc  ed 

□ 

Justlf loatlon - 

Bv  - — - 

Distribution/  _  _ 

Availability  Codes 
Avail  and/or 
Dist  1  Special 


/t 


UNCLASSIFIED 


MCUMITV  ekASMnCATieM  erTMW  PAWWItm  Vm* 


I.  INTRODUCTION 


» 

The  classical  work  of  Dodd  and  his  covorkers  at  the  Oak  Ridge  National 
Laboratory  deals  with  the  analysis,  design  and  optimization  of  eddy-current 
proble  colls  wound  around  an  air  core.  Hany  applications,  however,  require 
that  the  magnetic  field  produced  by  the  probe  coll  be  "shaped”  or  confined 
to  certain  regions  of  space,  especially  at  higher  frequencies,  and  this 
necessitates  the  use  of  highly  permeable  core  materials,  such  as  ferrites. 

4  The  problem  we  propose  to  solve  Is  depicted  In  Figure  1.  It  consists 

of  an  arbitrarily  shaped  body-o£-revolutlon  made  of  a  ferrite  (or  other 
highly  permeable  material),  excited  by  a  coaxial  coll,  and  in  the  presence 
of  a  plane-parallel  stratified  half-space  of  conducting  materials.  The 
half-space  simulates  the  workpiece,  in  eddy-current  nondestructive  evalua¬ 
tion  (NDE)  parlance.  As  Figure  1  Indicates,  the  system,  composed  of  the 
core,  coll  and  half-space.  Is  axlsymmetrlc;  this  is  the  simplest,  realistic, 
three-dimensional  problem  of  any  degree  of  generality. 

Though  the  model  that  we  will  develop  is  applicable  to  arbitrarily 
shaped  bodles-of-revolutlon,  we  will  specifically  attack  the  problems 
posed  by  the  core  shapes  of  Figures  2  and  3.  The  E-shaped  core  of  Figure  2 
is  typically  used  In  low-frequency  eddy-current  NDE  work,  whereas  the 
truncated  cylinder  of  Figure  3  Is  used  in  high-frequency  work.  Again, 
the  problems  posed  by  Figures  2  and  3  possess  an  axis  of  rotational 
synmetry. 

If  the  excitation  of  the  ferrite  is  weak  then  the  problem  may  be 
created  as  linear,  which  means,  among  ocher  things,  that  there  will  be 
no  harmonic  distortion  of  single-frequency  sinusoidal  signals.  At  the 
outset  of  this  investigation  we  make  the  small-slpial,  linearized 
assumption,  but  we  Intend  to  lay  out  the  appropriate  paths  to  follow  in 
^  dealing  with  the  more  general  nonlinear  problem. 

Cores  that  are  made  from  powdered  ferrite  materials  are  isotropic, 
but  those  made  from  single  crystal  ferrites  are  anisotropic.  The  math¬ 
ematical  approach  chat  we  use  in  developing  the  model  will  be  sufficiently 
general  that  it  could  Include  anisotropic  core  materials,  as  well. 

The  problem  solved  by  Dodd,  et  al.  [1-3]  differs  from  ours  in  that 
the  ferrite  core  is  absent.  Thus,  in  the  problem  they  solved  the 
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excitation  source  was  also  parallel  to  the  workpiece  (l.e.,  there  was 
no  Irregularly  shaped  core  to  destroy  the  plane-parallel  symmetry)  and 
this  permitted  them  to  solve  analytically  a  number  of  important  problems. 
We  can  follow  their  analytical  approach  only  so  far  (that  point  being 
the  determination  of  a  Green’s  function),  and  then  we  must  continue  with 
a  numerical  approach  In  order  to  compute  the  magnetization  within  the 
ferrite  core.  Having  confuted  the  core  magnetization,  we  can  compute 
the  fields  within  the  workpiece  analytically,  as  well  as  the  driving- 
point  impedance  of  the  coll/core  combination. 

In  a  sense,  Dodd  et  al.'s  work  brought  eddy-current  NDE  "of  age"  by 
showing  how  one  could  get  useful  analytical  results  bas.2d  on  a  rigorous 
application  of  electromagnetic  field  theory.  In  this  paper  we  extend  the 
application  of  electromagnetic  field  theory  to  Include  numerical  tech¬ 
niques  for  analyzing  the  fields  produced  by  complex  structures.  Our 
approach  Is  very  much  In  the  spirit  of  contemporary  research  In  problem¬ 
solving  In  electromagnetics  [4,5]. 

Ve  start  by  replacing  the  ferrite  core  by  an  equivalent  controlled 
source  of  Amperlan  currents,  which,  together  with  the  true  current  In 
the  exciting  coll,  coiiq>rlse  the  total  source  of  the  electromagnetic  field. 
The  field  is  expressed  as  an  Integral  over  the  regions  occupied  by  the 
source  currents,  l.e.,  the  core  and  coll.  The  integrand  is  a  vector 
function  of  two  arguments,  one  the  source  point  occupied  by  the  currents, 
and  the  other  argument  being  the  field  point  at  which  the  electromagnetic 
field  Is  to  be  evaluated.  The  Green’s  function  confuted  by  Dodd,  et  al., 
makes  Its  appearance  In  this  integral  expression. 

The  controlled  Amperlan  source  current  density  Is  not  known  a  priori, 
however,  becaiise  It  depends  on  the  value  of  the  field  at  the  source  point. 
Hence,  we  end  up  with  an  equation  whose  unknown  appears  both  outside  and 
inside  an  Integral  operator — an  integral  equation  for  the  unknown  Amperlan 
source  current  density  (which  Is  directly  related  to  the  magnetization 
of  the  ferrite  core).  This  integral  equation  Is  reduced  to  an  algebraic 
system  by  the  method  of  moments  [6],  and  Is  then  solved  using  a  linear 
equation  solver.  Having  the  Amperlan  and  true  currents,  we  can  compute 
the  electromagnetic  field  at  any  point  of  space.  Including  within  the 


workpiece,  by  straightforward  integration.  The  derivation  of  the  integral 
equation  and  other  analytical  matters  are  the  subject  of  Section  II. 

The  transformation  of  the  Integral  equation  by  means  of  the  method  of 
moments  is  dealt  with  in  Section  III,  and  Section  IV  discusses  some  of  the 
important  numerical  matters,  such  as  computation  of  special  functions  and 
Infinite  Integrals,  that  are  crucial  to  the  model. 

Results  of  the  model  are  presented  In  Section  V.  In  particular,  we 
show  the  approximate  7  field  within  the  two  cores,  the  driving-point 
Inqpedance  as  a  function  of  frequency  and  the  Induced  eddy-currents  as  a 
function  of  radial  coordinate  and  depth  within  the  workplace. 

A  section  on  covnents  and  extensions  of  the  model  concludes  the 
paper . 


II.  ANALYSIS 

(a)  Derivation  of  the  Volume  Integral  Equation 

In  the  sinusoidal  steady-stati.  Maxwell's  equations  are: 


V  X E  -  -jwB  (1) (a) 

V  X  H  -  J  (b) 

But  H  ■  E/Uq  -  M,  so  that  (l)(b)  becomes 

V xB  -  JuyQeQE  +  WqVxM  +  u^J  (2) 

In  (1)  and  (2)  E  Is  the  electric  field  Intensity,  B  the  ma^etlc 
Induction,  H  the  magnetic  field  Intensity,  M  the  magnetisation  (or  magnetic 
dipole-moment  per  unit  volume)  and  J  Is  the  true  electric  current  density. 
The  curl  of  M  in  (2)  clearly  plays  the  same  role  as  J  in  serving  as  a 
source  for  B;  this  source  term  is  the  Amperlan  current  density.  It  will 
play  a  fundamental  role  In  the  Integral  equation. 

Mow  let  B  ■  V  X  A,  E  *  -j(i£,  where  A  ■  Aa^  Is  the  vector  potential. 

It  hss  only  a  ^component  because  the  current  density,  J,  does  also.  Thus, 


A  (and  E)  are  aolanoldal  as  long  as  the  current  density  Is.  When  these 
new  definitions  are  substituted  Into  (2),  we  obtain  the  fundaaental  vector 
wave  equation  for  A: 


V  X  7  *  A  -  k^A  +  Mq  7  X  M  +  MgJ  (3) 

The  parameter  In  air;  within  the  workpiece  the  vector 

wave  equation  for  A  Is  the  same  as  (3) ,  except  that  H  >  0  and  k?  ■ 

2  ^ 
b)  V1q£q  ~  where  o  Is  the  electrical  conductivity  of  the  workpiece 

(we  are  assuming  that  the  workpiece  Is  a  homogeneous  half-space) .  The 
magnetic  permeability  of  the  workpiece  is  the  same  as  that  of  free-space, 
Uq,  which  means  that  the  workpiece  Is  nonmagnetic  (though  this  Is  not  an 
Important  restriction) ,  and  the  dielectric  constant  of  the  workpiece 
will  simply  be  assumed  to  be  e^.  Thoughout  this  development  we  Include 
the  displacement  current  density,  though  It  Is  negligible  at  the  fre¬ 
quencies  of  Interest,  even  In  air. 

We  begin  the  analysis  by  stating  the  vector  Green's  Identity: 


•7x7xQ-Q»7x7x  PldV  - 

^  [W  X  7  X  P  -  P  X  7  X  Q]  •  a  dS  (4) 

s  n 

where  the  vector  fields,  P  and  Q,  are  arbitrary,  except  that  they  must 
be  well-behaved  within  V  and  on  S.  We  are  going  to  apply  this  Identity 
to  certain  Green's  functions  and  vector  fields,  the  vector  fields  being 
the  A-flelds,  which  satisfy 

7  X  7  X  Aj^  -  kgAj^  -  11^7  x  M  +  7  x  7  x  A^  -  k^A^  ■  0,  (5) 

In  regions  1  and  2,  respectively  (see  Figure  4).  The  Green's  functions 
satisfy 

(r'lZ')  In  1 

V  »  V  .  0^  -  (r.r)  in  1  («(.> 


I 


7  X  7  X  -  0. 


(r' ,2')  in  1 
(Xyz)  In  2 


(b) 


j; 


Clearly,  the  Green's  function  Is  the  response  to  a  circular  filament 
of  current  of  radius  r*  located  on  the  plane  z-z'.  The  notation  G^^ 

(r,z;  r',z')  denotes  that  the  field  point  (r,z)  lies  in  region  i,  and 
the  source  point  (r',z')  lies  in  region  j,  as  shown  in  Figure  4. 

Now  apply  the  vector  Green's  identity  to  (5)  and  (6),  letting  P  <■  G. 
andG-A.  The  volume  of  Integration,  V,  must  be  broken  into  two  parts, 
corresponding  to  the  two  regions  of  Figure  4.  The  reason  for  this  is 
that  the  materials  differ  from  region  to  region,  and  so  the  field  vectors 
will  experience  certain  discontinuities  at  the  surface  separating  the  two 
regions  (z*0,  in  Figure  4).  Hence,  over  region  1  we  have 

.  7  .  V  .  V  X  M  +  .  J 

X*7x'7xG  »  k^  *0  +  "a  •  A. 

\  *  11  0^1  ®11  Slip  (f>  \ 

The  volume  Integrand  over  region  1  becomes,  therefore: 

where  Aj^  «  A^a^ .  Now  we  Integrate  with  respect  to  the  field  variables 
(r,z,^),  and  get: 

Vu  •  '  *  “  -  VGY  - 


#  [Aj  X  7  X  Gu  -  Gjj  X  7  X  Ajl  •  X^GS. 


(7) 


The  surface  consists  of  the  z-axls,  the  Infinite  spherical  sector, 
S^,  and  the  surface  z  0,  as  shown  in  Figure  4.  The  integral  along  the 
z-axls  vanishes  because  there  is  zero  surface  area  there,  and  we  assume 
that  there  exist  no  singular  sources  concentrated  along  the  axis.  The 
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Integral  along  vanishes  because  the  fields  and  satisfy  a 
radiation  condition  at  Infinity,  which  means  that  there  are  outgoing 
waves,  only,  at  Infinity.  This  leaves  only  the  Integral  along  z  ■  0, 
which  can  be  written: 


~2ir|o[Aj^  X  V  X  -  G^^^  x  V  x  1  •  a^  rdr  . 


z“0 


Hence,  we  have 


A^(r',z')  -  IitUq  11  J(r,z)  •  G^j^(r,z;  r’,z’)rdrdz  + 


coll 


+  I  V  X  M(r,z)  •  Gj^j^(r,z;  r',z')rdrdz  +  2ti* 


core 


'OV'- 


z)  xVxG^^(r,z;  r’,z*)  - 


-  G-,(r,z;  r',z')  x  V  x  A-(r,z)l  1  ‘a  rdr 

^  z-0  * 


(8) 


for  Che  result  of  the  integration  over  V^,  assuming  that  (r',z')  lies  In 
(recall  the  meaning  of  G^^ :  both  field  and  source  points  lie  In 
Note  that  both  the  core  and  coll  lie  In  region  1. 

Now  we  continue  with  the  application  of  Green's  vector  Identity  to 
region  2.  Here  we  use  the  second  equation  In  (5),  together  with  (6)(b). 
Keeping  In  mind  Chat  the  source  point,  (r’,z'),  still  lies  within  region  1, 
we  get,  after  arguing  as  before: 


0  -  2Tr|* lA2(r,z)  X  7  X  G2j^(r,z;  r',z')  -  G2j^(r,z;  r’,z')  x  7  x  A^Cr.z] 


z»u 


•  a  rdr 
z 


(9) 


This  integral,  however.  Is  proportional  to  the  third  integral  In  (8), 
as  the  following  argument  shews.  Use  the  scalar  triple-product  to  rewrite 
the  Integrand  In  (9)  as 


(10) 


I(a^  X  k^)  •  V  X  -  (a^  x  G2j^)  •  V  x  A^) 


Now  a  X  A  Is  the  tangential  component  of  A,  which  is  proportional 
z  _ 

to  the  tangential  coiiq>onent  of  E,  which  must  be  continuous  at  the  boundry, 
z  ■  0.  Hence, 


a  X  A-  *  a  .  X  A, 

z  2  z  1 


at  z  *  0.  This  is  also  the  first  boundary  condition  to  be  Imposed  on  G: 


a  X  G--  -  a  X  G,, 
z  21  z  11 


at  z  0.  Therefore,  (10)  becomes 


[(a  X  L)  •  7  X  G  .  -  (a  x  G  )  •  V  x  A  1  1 
z  1  21  2  11  2  ^,0 


[-(a^  X  7  X  G^j^)  •  A^  +  (a^  X  7  X  A^)  *  I  » 

2«0 


where  we  have  again  invoked  the  scalar  triple-product.  We  note  that 


■—  a  X  7  X  A  1  -  ~  a  X  7  X  A.  1  , 

*  2^.0  Vq  z 


which  expresses  the  continuity  of  the  tangential  H-flelds  at  z  ■  0.  Again, 
we  impose  the  same  boundary  condition  on  6: 


^  a  X  7  X  ?  1  -  ^  **  »  ^  *  ^11  I 

*^1  *  ^^z-0  *  ^^z-0 


and  when  this  is  put  into  the  right-hand  aide  of  (11) ,  we  get 


i 


7 


‘-<S  X  V  X  G^)  •  Xj  +  X  7  X  A^)  .  5^1 


i  [A^  X  7  X  0^  -  0^  X  7  X  A^l  jj  5^ 


which  Is  times  the  lotegrand  of  the  third  Integral  in  (8) ,  as  we 

set  out  to  prove.  Hence,  because  (9)  vanishes,  so  does  the  third  Integral 
in  (8). 

To  susmarlze,  therefore,  upon  imposing  the  boundary  conditions  at 
*  -  0: 


*  SiCr.z;  r’.*')  -  5^  X  G^^(r,z;  r*,z') 


(22)1 


^  S  X  V  X  G,,(r,z;  r’,z’)  -  x  V  x  G,,(r,z;  r’,z’) 


^  * 


^0  * 


then  (8)  becomes 


z’)  -  2'ityQ  j|  J(t,z)  •  Gj^j^(r,z;  r',z')rdrdz 


+  2itVq  II  V  X  M(r,z)  •  Gj^j^(r,z;  r',z')rdrdz  . 


We  eaq>haalze  that  the  differential  operation  in  (12) (b)  is  with 
respect  to  the  field  variables,  (r,z).  From  here  on  we  take 
as  we  stated  earlier. 

The  magnetic  Induction  is  the  curl  of  A: 
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^  B(r,z)  -  7  X  A(r,2) 

-  V  X  II  Gj^j^(r,z;  r'z')  •  j(r' ,z')r’dr'dz'| 

coll 

+  V  X  |a^2TnjQ  Gj^j^(r,z;  r',z')»7’  x  M(r' ,z’)r’dr'd2'|-  ,  (14) 

core 

^  where  we  have  used  the  symetry  of  G^^: 

GjjCt.zi  r'.z')  -  OjiCt'.z':  r.z)  . 

y  and  we  note  that  the  curl  operation  is  with  respect  to  the  primed  variables 

under  the  Integral  sign,  and  with  respect  to  the  field  variables  (the 

unprined  variables)  outside  the  integral  sign. 

The  first  term  of  (14)  is  the  incident  field,  due  to  the  coil, 

-(s) 

and  the  second  is  the  acattered  field,  B  ,  due  to  the  induced  polarization, 
-  '•(s) 

M.  We  proceed  to  reduce  B  .  Start  with  the  vector  identity 

Gj^j^(r,z;  r',z')  •  7*  xM(r',z’)  -  -7'  •  r',z’)  xM(r',z')j- 

+  M(r',z')  •  7’  X  Gj^(r,z;  r’,z’)  ,  (15) 

and  write  B  (r,z)  as: 

« 

B^®^(r,z)  •  2'rru^7x  |a^  ||  M(r',z')  •  7'  x  Gj^j^(r,2;  r' ,z')rMr'dz' 

core 

f  II  7’  •  lG^^(r,z;  r'.z')  x  M(r’ ,z’)]rMr’dz’J.  (16) 

core 

Upon  using  Gauss'  theorem  on  the  second  integral,  and  extending  the 
^  surface  of  the  core  slightly  beyond  the  (finite)  support  of  M  (H  is 

smooth,  so  that  Gauss'  theorem  is  applicable),  the  second  integral  van¬ 
ishes.  We  are  left  with 
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B^®^(r,z)  -  2irviQVx 


II  M(r',z')  •  V  X  Gj^^(r,z;  r' ,z')r*dr'dz'  (17) 


core 


The  total  field  satisfies: 


B(r,z)  -  B^^^(r,z)  +  B^®^(r,z) 


-  B^^^(r,z)  +  2iruQ  V  X  a^  II  M(r’,z’)  •  7  x  G^j^(r,z;  r'z')r'dr’dz 

(18) 


core 


To  make  this  into  a  volume  Integral  equation  for  B  or  H,  we  note  that 


fi  .  5/v„  -  S  -  i(^  -  i)  , 


where  )i  is  the  permeability  of  the  ferrite,  so  that  (18)  yields  two  equiva¬ 
lent  integral  equations: 

B(r,z)  -  B^^^(r,z)  +  2ir  7  x  ||  (1-  •^)  B(r',z’)  •  7'  x 

core 


X  Gj^j^(r,z' ;r' ,z')r'dr'dz' 


(19) (a) 


(•jj;~)M(r,z)  •  B^^^(r,z)  +  27ryQ  7  x  a^  || 

®  core 


M(r*,z')  •  7'  X  Gj^j^(r,z;  r*  ,z')r'dr*dz’ . 


(b) 


In  either  equation  the  unknown  vector  field  can  be  expanded  in 
vector  pulse-functions,  because  the  curl  operation  has  been  transferred 
tc  the  Green's  function. 

The  curl  operation  In  the  unprlmed  coordinates  will  be  "mollified," 
or  smoothed,  by  using  an  appropriate  test-function  in  the  method  of 
moments.  To  see  how  this  is  accomplished,  consider,  again,  the  vector 
identity  that  was  used  in  (15) : 


10 


P«7xQ--V-(PxQ)+Q*7xP  . 


(20) 


If  P  satisfies  V  x  P  -  0,  then  (20)  becomes  P  •  V  x  Q  -  -V  »  (P  x  Q) . 
When  this  result  is  integrated  over  a  patch  of  space,  and  Gauss'  theorem 
invoked,  we  get 


II  [P  X  Q1  •  5^dS 
dPatch 


(21) 


This  result  is  void  of  differential  operators.  This  is  i^t  is 
meant  by  "mollifying"  an  equation.  We  will  use  this  on  (19) ,  treating 
the  bracketed  term  on  the  right-hand  side  as  the  vector  Q. 


(b)  Derivation  of  the  Green's  Function 

The  Green's  function  that  appears  in  the  integral  equations,  (19), 
satisfies  (6)  and  (12) .  Because  there  is  only  a  (|>-coinponent  of  6,  these 
equations,  in  cylindrical  coordinates,  reduce  to: 


6(r-r')6(z-z' 


-,  0<z,z'  (22)  (a) 


”'°21  .  ”'”21  .  1 
3,1  3,1  "  1' 


21  2 

r  ''l°21  “  ° 


z  <  0 


O^Cr.O,  r'..')  .  Oj^Cr.O,  r'.r-) 


(r,0:  t',,')  . 


''*21 

(r,0;  r',z') 


In  going  from  (12) (b)  to  (22) (d)  we  have  put  •  y^. 

In  solving  (22) (a)  we  must  consider  two  subregions  of  region  1: 
z  >  z',  z  <  z'.  At  z  •  z'  we  require  that  the  Green's  function  be 
continuous,  but  that  there  be  a  discontinuity  in  its  slope.  The  amount 
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of  discontinuity  can  be  determined  by  Integrating  (22) (a)  over  a  small 
Interval  about  z ' : 


fz'+e  1 

[Left-hand  aide  of  (22)(a)]dz  =  | - ^  |  - - 


Thtis,  the  discontinuity  in  slope  is  given  by 


3z  'z| 


The  solution  of  (22)  Is  facilitated  by  using  Bessel  transforms.  We 
start  by  stating  the  fact  that  any  function,  f(r),  can  be  expanded  as 
the  double  Integral 


(rDldZ 


rf(p)j^ 


(£p)pdp 


f^lch  means  that 


-  Tj  (r^)j  (plyldl 
P  io  “  “ 


The  second  Integral  In  (24)  defines  the  Bessel  transform,  F(4), 
of  f(r). 

We  apply  the  Bessel  transform  to  (22) (a)  by  writing 


Gii(r,z;  r',z')  -  j5j^j^(£,z;  r' ,z')J^(^r)^ 


+  -  feiiU,*;  r\z’)l^iy'(lr)’Kr^CZt)/Zr 

3'  '  n 


-  j(ix)/arriiiu 

a 
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But  J^(z)  satisfies  Bessel's  equation 


+  j;_(*)/z  -  J^(z)/z^  -  -Jj^(z)  . 


Hence,  in  the  last  Integral,  and  hereafter,  ve  let  n  -  1  and  get 


!u 


^  Giia,z;  r’,z•)JJ^(^r)^d^ 


The  z-derlvatlves  are  simply 


r 

- (£,z;  r',z')J,(£r)£de  . 

32  Jq  dz  ^ 


Upon  assembling  these  results,  together  with  (25),  Into  (22) (a),  we 


2^  ^00 

fflin  +  j  J^(r^)J^(r’^)^dC  .  (27) 

Jq  dz^  ° 


Since  this  holds  for  all  r  we  may  simply  equate  Integrands  to  get 
an  ordinary  differential  equation  for  the  Bessel-transformed  Green's 
function: 


d\  2  2 

- -  -  «(z-z’)  ,  0  <  z,z'  . 

dz 


(28) (a) 


Similarly 


^  +  ik\-Lh^21  “  °  * 


z  <  0,  z'  >  0  . 
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The  boundary  condltlona  (22)(c),(d)  and  the  discontinuity  condition, 
(23),  when  expressed  In  terns  of  the  Bessel  transforms  are 


(c) 


d6.-  dCj- 

(£.0;  r',z')  -  (£,0;  r',z') 


dz 


(d) 


d6  1  dS. .  —  J 

r'.z')  ^  (^.*1;  r',z')  -  (r^)  . 


(e) 


Also,  we  have  continuity  of  at  z  z* : 


r'.z')  .SyU..',  r'.r') 


(f) 


Solutions  of  (28) (a)  in  the  appropriate  subregiuns  are: 


-a^(*“z') 

6j^j^(^,z;  r*,z')  •  Ae  ® 


,  z  >  z'  >  0 


-a  (z-z')  a-(z-z') 

-  Be  +  Ce  ,  z'  >  z  >  0  (29) 

2  2  1/2 

where  >  (£  -k^)  .  Equations  (28) (e) , (f )  are  applied  to  determine 

relations  between  the  arbitrary  constants  A,  B,  C: 


A  -  B  -  C  -  0 

J,(r'^) 

^  ®  ^ - 'iisr-' 


(30)  (a) 
(b) 


The  solution  of  (28) (b)  requires  only  a  single  arbitrary  constant: 


a  (z-z') 

e2^(^.z;  r',z')  -  De  ^ 


z  <  0,  0  <  z' 


(31) 
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where  When  (28)(c),(d)  are  applied  to  (29)  and  (31), 

we  get  the  final  two  equations 


— ot  z*  ot^z'  ~ol  z* 
1  0  0 
De-^  -Be  +Ce“ 


(32) (a) 


a.r'  -a  z' 

-a^Be  +  a^Ce  ^  -  a^^De  ^ 


(b) 


The  solution  of  (30)  and  (32)  Is  readily  found  to  be 


01^  —  a  “20  z’ 

0  “l  ^  0  ^  ^ 


47ra^ 


(33) (a) 


“o““i 


-2“0*’ 


Airci, 


(b) 


Jj^(r’£) 

4TO„ 


(c) 


-(a  -a  )z’ 
J^(r’£)e  0  ^ 

°  “  (aJj5o!]72i 


(d) 


The  final  expression  for  the  Green's  functions  Is  obtained  by 
substituting  (33)  Into  (29)  and  (31) ,  and  then  that  result  Into  (26) : 


®u<' 


..i  i 

n  L 


-a-(z+z')  -a.(z-z'J 

Re  +  e 


2a, 


JJ^(r'^)JJ^(r^)^d^,  z>*'>0 


-o„(t+z')  aQ(z-z')-i 
e 


-  ; 

1  r  S£ - +j 

-  ^  J  [  2a, 


Jj^(r'^)Jj^(r-£)'td^»  z'>s>0 
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I 


<»r  -a  (z+*‘)  -a  |«-z' 

Re  +  e  “ 


J^(T*l)J^(Tl)Zdl,  z,z*>0  (34) (a) 


(r.z;  r'.z')  ^ 


(a.z-a.z') 


J^(r’i)J^(r£)£d£ 


We  have  defined 


p  I  ■■  ^  T  «  ■  ■  ■ 

“o-^r  “o-^i 


In  (34).  We  have  already  alluded  to  the  fact  that  6^^  Is  symmetrical  under 
Interchange  of  source  and  field  coordinates,  and  (34) (a)  verifies  that 
fact.  The  final  form  of  6^^  in  (34) (a)  is  merely  a  compact  representation 
that  Includes  the  first  two. 

Now  that  we  have  the  Green's  function  that  is  to  be  used  in  (19), 
we  can  proceed  to  transform  the  Integral  equation  into  vector-matrix  form. 


Ill.  REDUCTION  OF  THE  INTEGRAL  EQUATION  TO 
VECTOR-MATRIX  FORM:  THE  METHOD  OF  MOMENTS 

From  here  on  we  will  work  with  the  Integral  equation  (19) (a).  The 
starting  point  in  the  numerical  reduction  of  this  equation  is  to  parti¬ 
tion  the  core  In  a  regular  grid,  as  shown  In  Figure  5  for  the  cylinder. 
Based  on  this  partition,  we  expand  the  unknown  field  In  pulse  functions 

»c_ 

B(r,*)  -  E  B.P^(r,*),  (3 


or  In  component  form 
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(37)(*) 


J-1  J  'b  J 


B  (r.z)  -  E  bJ*^P.(r,*)  . 

*  4_1  J  J 


where  {b^*^},  {b^*^}  ere  the  expansion  coefficients  for  the  r-  and  z- 
components,  Is  the  number  of  cells  In  the  grid  and  P^(r,z)  Is  the  jth 
pulse  function,  which  Is  defined  by 


Pj(r,z)  -  1,  (r,z)  in  jth  cell 

■  0,  otherwise. 


Upon  substituting  (36)  into  (19) (a) ,  we  get 


Nc 

E  B.P.(r,z)  •  B^^\r,z)  +  E  7  x  a.27r  [[ 

1-1  J  J  1-1  ’ 

J  J  jth  cell 


(1- y)Bj(r’,z’)  «  7’  xG^^(r,z;  r' ,z')r’dr'dz' .  (39) 


The  next  step  In  the  method  of  moments  is  to  take  moments  of  (39) , 

l.e.,  to  multiply  (39)  by  testing-functions  and  then  Integrate.  For  our 

testing-functions  we  take  the  sazM  functions,  VP^(r,z),  where  V  Is 

either  (r/r.  )a  or  a  ,  as  for  expansion,  (36).  Wien  the  same  functions 
D  r  z 

are  used  for  expansion  and  testing,  the  awthod  Is  called  the  Galerkln 
variant  of  the  method  of  moments.  Thus,  the  result  of  taking  mmaents  of 
(39)  Is 


E  II  •  TjPj  (r,z)rdrdz  • 


V  •  B^^^(r,z)  rdrdz  + 


ith  cell 


1th  cell 


E  j|  V  •  7  x{a^fj(r,z)}rdrdz, 
ith  cell 


t 
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where  we  have  written 


» 


I 


t 


t 


f^(r,*)  -  2it  JI  (l~)B^(r',*’)  •  7'  x  Gj^^^Cr.z;  r' ,z’)r’dr'd*' 
jth  cell 


■  ff 


-3G^^(r,z;  r',z') 


3z' 


-)r'dr'dz' 


jth  cell 


ff  ^0  r^^ll^*'***  6n(r.2;  r’.z’)] 

II  — j' 


'dr'dz'. 


jth  cell 


Recalling  the  definition  of  allows  us  to  reduce  the  left-hand  side 


of  (AO)  to 


‘ij  II  ’  •  *1’'“*'“ 


(42) 


ith  cell 


The  first  Integral  on  the  right-hand  side  of  (40)  becomes  the  forcing- 
function  vector,  fj. 

Because  V  x  V  -  0,  the  second  integral  on  the  right-hand  side  of  (40) 
may  be  reduced  by  the  use  of  (21),  with  the  result  that  the  Integral 
becomes 

I  •  a^fj(r,z)dA  , 

ith  cell 

where  the  contour  integral  ii^lies  an  integral  over  the  boundary  of  a 
two-dimensional  cell.  Upon  referring  to  Figure  6,  this  integral  may  be 
written 


f 
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(41) 


9 


^  [l^  X  V]  •  S^f^(r,z)dA-J  )(a^  * 


1th  cell 


*1  f'i 

+ 1 , /j  ‘’'r'***  *  ■  V*  *  ’’ 

*1  'i 

■  1  *  V*- 


Now  we  let  V  •  (r/r.  )a  in  (43),  and  get 

D  t 


I  [a^  X  i,]  •  a^(’^/»‘b>f 


1th  cell 


f'l 

(r,z)dA  - 

^  i(-) 


(r/r^)f^(r,2^“^)rdr 


•  / 

-f  (r/r  )f  (r,8{'*’^)rdr 

J  (-)  ^  ^  ^ 


Similarly,  when  V  *  a 


t  '*♦  *  *«’ 


1th  cell 


*1  ^ 


z) 

l(-> 


Upon  eubetltutlng  (41)  into  (44) (a), (h),  we  get,  reepectiwely 


1 


We  can  now  start  to  assaable  the  matrix  equations.  When  V  -  (r/r.  )ar 

l> 

in  C^O)  we  get,  upon  using  (42)  and  (45): 


.(+) 


"i'*  ll  II  • 

ith  cell  ° 

•  r'.z’)  -  G^(r,^^■^^  r',z')]  J(r’ ,z')r'dr'dz' 

N 

+  z'bj'^27r  [  (^)rdr  [f  (i--^)(|l)  . 

^  iJLii  “  '» 


3z' 


3z' 


r’dr'dz* 


.(+) 


+  E°bJ''2itf*  (^)rdr  f(  (l-^)- 
i<->  ■>  Jtb'i.11 

r30^(r,z'->j  t'.z')  Ojj(r.z<->:  r'.z') 

•|_ - p - + - p - 


3G 
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(r.zj'*’^;  r’.z')  G  (r.zj'*’^;  r',z’)l 
- 


dz* 

1-1,  ...,  N. 


(46) (a) 


and  when  V  ■  a  ,  (40)  becomes 

X 
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r'.z')  -  r'.z')!  J(r' ,z')r'dr'dz' 


H  ^(+) 

If  <^-T>t>- 

J  '  4t-U  ° 


2^  '  jth  cell 


(_)  '*»*’)  (+)  r’.z') 

■1  a?  'i  a? 


r’dr'dz' 


N  z^"^^ 

f  '  dz  ff 

J  •*•  4#.f.  A. 


(l-f) 


Jth  cell 


fr . 


3Gii(r<+>.z;  r*,z')  G^^(r[-^>.z;  t'.z’) 

+ - - 


"sp- 


-  r 


(-) 


^^ll^'i'^’*’  G,,(r/"\z;  r',z') 


IP 


11' 


'  '  ,z;  r- ,z  )  \ 
i - p -  I  r'dr'dz' 


i  ■  1,  . . . ,  N  ,  (46) (b) 

c 


where  the  first  term  on  the  right-hand  side  of  these  equations  comes  from 
substituting  In  (14),  Into  the  first  term  on  the  r Ight-hand  side  of 

(40),  and  then  using  (21)  to  mollify  the  term,  as  before. 

These  equations  can  be  put  into  the  following  matrix  form: 


22 


I 


A^rr)  .  Q<irx) 


-  G 


(zr) 


,(zz) 


' 

1 

1 

N  1 
U  \ 
^  1 

~(r) 

-  G<“>J 

Jz) 

(47) 


The  overall  system-matrix  Is  of  order  x  2N^. 

The  elements  of  the  submatrices  and  forcing-function  vectors  can  be 
reduced  to  infinite  Integrals  by  substituting  the  Green's  function, 

(34) (a),  and  its  curl-components,  into  (46).  The  components  of  ^  G^ 
are  given  by 
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3Gii(r,2;  r'.z') 
3z* 


_1_ 

2tt 


r 

0 


-Re 


3r'  r* 


K  -“qIz-zT 

+  e 

2 

J 

(-) :  z'  >  z 

(+):  z  >  z' 

_L  r 

2S  J 

0 

2«0  J 

(48) (a) 


(b) 


In  evaluating  the  Integrals  that  define  the  matrix  elements,  we 
interchange  the  orders  of  Integration  with  respect  to  (r,z,  r',z')  in 
(46)  and  those  with  respect  to  the  transform  variable,  t.  The  results  are: 


,(r).W 


0 


2Fi(*i»VV*i  ’*i  25: 


dl 


(49) (a) 
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■  '  '** 


(i)^dr  .  dyrX-)-z<  >) £ 

.(-)  -(-) 


.(+)  ,<+) 


Ij 


dz 


.(-)  .(') 


I  ,d, .  5,  .re<z«-z<->)  •  'M .:  fli  -5L 


'lj‘b'“l  -1 


The  functions  used  In  (49)  are  defined  here: 


c 

?j^(a,b)  •»  I 


-a^a  -a„b 
e  "  e  0 

e  dz  ■  - r - 


(50) (a) 


F^Ca.b.c)  ■  j  e 


-a„|c-.|  .-“o'-')  .  .-“o'»-') 

dz  - - — ~”r - 


c  <  a 


-a  (c-a)  -o  (b-c) 
2-e  ^  -e  ^ 


a  <  c  <  b 


-a-(c-b)  -o-(c-a) 


T 


t 


-a  (c-a)  -a  (b-c) 
e  -  e 


a  <  c  <  b 


-OgCc-a)  -a  (c-b) 
e  -  e 


b  <  c 


(d) 


b  d  ^  I  ii 

-ajz-z’' 


FgCa.b.c.d)  -  dz  1  dz'  (+  e  “  )  , 


a  c 


(+) ;  z  <  z' 
(-) :  z'  <  z 


F^Ca.b.c)  -  F^Ca.b.d) 


ct. 


(e) 


Bj^Ca.b) 


zjQ(z)dz  -  bJ^(b)  -  aJj^(a) 


(f) 


D 

B^Ca.b)  -  I  z^J^(z)dz  -  b^J^Cb)  -  a^J^Ca) 


(g) 


I  (a 


t 

,b)  -  j 


zJj^(z)dz 


(h) 


We  will  discuss  the  computation  of  the  last  three  functions  In  the 
next  section. 

An  Interesting  feature  of  the  matrices  In  (49)  Is  their  sTimiietry 
relations: 


p(rr)  (rr)  -(zz)  -(zz)  „(rz)  _(zr) 

®1J  *  ®J1  '  ®1J  ®jl  ’  ®1J  ®jl 


(51) 


These  relations,  which  reduce  the  time  to  compute  and  fill  the 
system-matrix,  are  a  result  of  uslx:g  the  Galerkln  variant  of  the  method 
of  moments. 
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In  deriving  the  forcing-function  vector  in  (49) (a) , (b) ,  we  assumed 
that  the  coil  is  uniformly  and  densely  wound  with  n^  turns  per  unit  area, 
that  it  occupies  the  rectangle  r^  <  r  <  r^,  <  z  <  Z2,  and  that  it 

carries  a  current  of  I  amperes. 


IV.  NUHERICAL  ANALYSIS:  COMPUTATION  OF 
SPECIAL  FUNCTIONS  AND  INTEGRALS 

The  method  of  moments,  or  any  other  scheme  that  is  used  to  generate 
numbers,  depends  critically  on  effective  numerical  methods.  The  end  result 
of  the  method  of  moments  is  the  vector-matrix  equation  (47) .  We  solve 
this  equation  very  efficiently  by  using  a  linear  equation  solver  that 
employs  the  Crout  algorithm  to  factor  the  system-matrix  (we  do  not  use 
matrix  inversion,  because  that  is  less  efficient) .  The  algorithm  uses 
scaled  partial  pivoting  and  back  substitution;  it  is  implemented  in  complex 
double-precision  on  the  computer. 


(a)  Computation  of  Special  Functions 

The  Bessel  functions,  Jq(z)>  appear  throughout 

the  vector-matrix  elements  of  (49) ,  are  easily  computed  by  appealing  to 
a  series  expansion  when  abs(z)  <  10,  and  to  well  known  asymptotic 
expansions  when  abs(z)  >  10  [7,  Ch.  9]. 

The  computation  of  the  Bessel-function  Integral,  (50) (h) ,  is  based 
on  the  fact  that  [8,  Vol.  1,  p.  219]: 


I  t“j^(t)dt 
0 


mtn+l 

z _ 

2“(mfn+l)  r  (ni-1) 


F 

12 


(m+n+1)  /  2 
(m+ii+3)/2,  ttfl 


(52) 


where  z  is  a  complex  number  and  ^F^  is  a  generalized  hypergeometrlc 
function.  Here,  m  ■  n  ■  1,  so  that 


3/2 

2 

J  tJ^(t)dt-^  ,F2 

n 

5/2.  2 

-z 

4 

, 

(53) 
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t 


I 


Now 


f  3/2 


1^2 


5/2, 


2 


(54) 


vhere 


(2)j^  k  I  (2k+3) 


(55) 


Hence , 


Vn  (k+2)(k+l)(k  +  5/2)  *k 

is  the  recursion  relation  to  be  used  for  the  series  expansion,  (54),  when 
abs(z)  <  10. 

For  aba(z)  >  10,  we  need  an  asymptotic  expansion  for  (54).  Again, 
from  [8,  p.  199]  we  have 


1^2 


' 

3/2 

2 

-z 

5/2,  2 

4 

4 

(2ir)^/^(z)^^^ 


.  I  ,j(*-5iT/4)  ,  j.  ^_jjk  ^^^k  ^  ^-j(*-5Tr/4)  ^  (57) 

where  the  coefficients,  satisfy  the  recursion  relation 

djj  -  1,  d^  -  -7/8,  2(k-*-l)dj^j  -  (3k^  +  2k  -  7/4)d^  - 

-  (k  +  3/2)  (k  -  1/2)  (k  -  3/2)dj^_j^  .  (58) 

The  final  fora  of  the  asymptotic  expansion  Is,  therefore: 
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f  '  d./z'‘  + 


i^/z^  }  • 


(59) 


The  additional  term  of  unity  comes  from  the  fact  that 


«■  * 
tJ^(t)dt  -  -  zJq(z)  +1  Jq 


(t)dt 


-  -  zJq(z)  -  I 


jQ(t)dt  + 


JQ(t)dt 


-  zJqCz)  -  I 


JjjCOdt 


(60) 


which  agrees  with  (59). 

(b)  Numerical  Evaluation  of  Infinite  Integrals 

The  computation  of  the  Infinite  Integrals  that  appear  In  the  matrix 
elements,  (49),  la  performed  by  computing  a  sequence  of  partial  sums  of 
Integrals  over  finite  Intervals  and  then  applying  the  e-algorlthm  In  an 
attempt  to  accelerate  the  convergence  of  the  sequence  [9].  The  finite 
Integrals  are  computed  by  using  a  30-polnt  Gaussian  quadrature;  the  length 
of  the  finite  Interval  Is  determined  empirically.  The  method  has  been 
used  In  other  nuamrlcal  electromagnetic  codes  [10]. 

V.  RESULTS 

After  (47)  Is  solved  for  the  expansion  coefficients,  {b^'^}, 

Che  results  are  then  substituted  Into  (37)  to  produce  the  components  of 
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the  Induction  field,  B,  within  the  core.  With  this  field  in  hand  it  is 
possible  to  compute  other  functions.  In  this  section  we  Illustrate  some 
of  the  results  that  have  been  obtained  for  the  E-shaped  and  cylindrical 
cores.  Some  of  the  results  correspond  to  the  workpiece  being  absent,  in 
which  case  and  are  real.  When  the  workpiece  is  present 

{  b  and  {b^*^}  are  complex. 

(a)  Approximate  B-fleld  Within  the  Core 

Figure  7  shov/s  the  approximate  B-field  within  the  cylindrical  core, 
and  Figure  8  shows  the  same  thing  for  the  E>shaped  core,  both  when  the 
workpiece  is  absent.  In  each  case  the  line  segments  represent  the  approxi¬ 
mate  field  at  the  center  of  each  cell  of  the  grid  that  was  used  in  the 
method  of  moments.  The  line  segments  are  shown  emanating  from  the  bottom 
or  top  line  of  a  cell  simply  to  facilitate  the  presentation. 

The  radius  of  the  cylinder  is  0.0625"  and  the  height  is  0.500". 

The  inner  radius  of  the  coll  is  also  0.0625",  the  outer  radius,  0.125", 
and  the  height  is  0.250".  The  outer  radius  of  the  E-shaped  core  is 
0.1875",  and  the  height  is  0.500".  The  same  coll  is  used  with  the  E- 
shaped  core  as  with  the  cylinder.  The  permeablity  of  the  core  is  1000 

We  note  that  B  is  generally  smaller  in  magnitude  for  larger  radii. 

This  means  that,  because  the  cross-sectional  areas  increase  with  in¬ 
creasing  radius,  the  net  magnetic  flux  entering  the  region  is  equal  to 
the  net  flux  leaving,  as  is  to  be  expected  from  the  condition  V  ■  B  «  0. 
Indeed,  we  have  computed  an  approximate  value  of  7  ■  B  at  one  point 
within  the  cylindrical  core  and  have  found  it  to  be  zero,  certainly  within 
the  accuracy  of  the  expansion  in  pulse  functions.  This  is  a  very  desir¬ 
able  feature  of  our  results;  it  shows  that  the  model  produces  physically 
consistent  and  accurate  results  with  the  grids  shown. 

(b)  Driving-point  Impedance  of  the  Coil  and  Core 

The  driving-point  impedance  of  the  coll  is  defined  to  be  the  ratio 
of  the  coil's  induced  voltage  to  its  current.  The  voltage  Induced  into 
a  single  turn  of  the  coll  located  at  (r',z')  is  given  by 


V(t*,**)  -  -  I  f  •  dZ  -  jU)2irr'A(r',z’)  , 


(61) 


idvlch  follovs  from  the  fundameatal  relation  E  -  juA*  The  coll  is  uni¬ 
formly  and  densely  wound  with  n  turns  per  unit  area  and  occupies  the 

c 

rectangle  <.  r  £  ^2*  *i  5.  *  £  *2*  before.  Hence,  the  net  voltage 
Induced  Into  It  Is  given  by 


2  f  2 


V  ■  jcj2nti 


11 

'l  *1 


rA(r,z)drdz  . 


(62) 


The  coll,  of  course,  is  in  region  1,  and  the  vector  potential 
throughout  region  1  Is  given  by  (13),  which,  after  using  (15)  to  trans- 

f er  the  curl  to  the  Green' s  function  and  then  substituting  M  ■  B  (■; - -) , 

becomes 


A(r,z)  -  2TruQ  ||  r’ ,z’)  J(r' ,z')r'dr'dz' 

coil 


11 


+  2u  11  (1  --r)B(r',z')  •  r  X  G„(r,z;  r' ,z')r'dr'dz'  .  (63) 


11' 


cote 


Thus,  the  induced  voltage  Is  given  by  substituting  (63)  into  (62) 


*^2  *2  '2 


V  •  jw(2ir)  ^n^I^Q  |  rdr  |  dz  |  r'dt'  | 


dz'Gj^l(r,z;  r',z’) 


*1  *1  ’'l 


+  J<i)(2ir)^n^  I  rdr  |  dz  ||  (1  -  ~)B(r',z’) 


Zj  core 


7'  *  Gj^(r,z;  r',z’)r'dr'dz'  . 


(6A) 
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The  first  term  on  the  right-hand  side  of  (64)  Is  the  voltage  that  Is 
Induced  Into  the  coll  In  the  absence  of  the  core  (the  coll  self-voltage), 
and  the  second  term  Is  the  additional  voltage  Induced  Into  the  coll  due 
to  the  Induced  magnetization  within  the  core.  Upon  using  the  expression 
for  the  Green's  function,  (34) (a),  and  then  carrying  out  the  Integral  In 
the  first  term  of  (64) ,  the  coll  self-voltage  becomes 


5  RF^(z,  ,z-)  +  F  (z  ,z,,z  ,z,) 

J  - it  . 

0  L  0  J  ^ 


where  F^,  F^,  and  7  have  been  defined  In  (50). 

The  second  term  of  (64)  Is  evaluated  after  substituting  (48)  and 
(37),  and  then  evaluating  the  Integrals  with  respect  to  (r,z)  and  (r',z'). 
The  result  la 


Jw27r(l  - 


where  Fj^^Vy^n^I,  F^^Vugn^I  are  the  Integral  expressions  of  (49) (a), (b). 

The  total  voltage  Induced  Into  the  coll  is  the  sum  of  (65)  and  (66) , 
and  when  this  Is  divided  by  the  coll  current,  1,  the  result  Is  the  driving- 
point  impedance. 

We  have  computed  the  driving-point  l]q>edance  at  four  freuqencles, 
and  present  the  results  In  Table  1.  The  columns  In  the  table  correspond 
to  the  dlrvlng-polnt  Impedance  of  the  coll  alone,  coll'fcyllndrlcal  core, 
and  coU+E-shaped  core.  The  data  are  In  ohms,  and  correspond  to  the  dimen¬ 
sions  stated  In  Subsection  (a).  The  numbers  in  parentheses  are  powers  of 
ten  that  multiply  the  adjacent  number:  for  exaaiple,  j0.2958(-2)  ■ 

JO. 2958  X  10  As  a  benchmark  calculation  to  test  the  accuracy  of  our 
Integration  subroutine  that  uses  the  e-algorlthm  (recall  Section  IV. (b)), 
we  computed  the  Inductance  of  the  coll  In  air  and  compare  It  with  the 
classical  work  of  Grover  [11];  the  results  are  In  excellent  agreement. 

The  Inductance  of  the  coll  in  air  Is  determined  from  the  first  column 
In  the  upper-half  of  the  table. 


33 


TABLE  1:  DRIVING-POINT  IMPEDANCE  AT  FOUR  FREQUENCIES 


AIR 


COIL 

CYLXND&R 

ESHAPE 

FREQ 

j0.2958(-2) 

j0.1554(-l) 

j0.3271(-l) 

IkHz 

J0.2958(-l) 

j0.1554(-0) 

j0.3271(-0) 

lOkHz 

j0.2958(-0) 

J0.1554(+l) 

j0.3271(+l) 

lOOfcHz 

j0.2958(+l) 

J0.1554(+2) 

j0.3271(+2) 

IMHz 

WORKPIECE 

COIL  X(-3) 

CYLINDER  X(-3) 

ESHAPE  X(-3) 

FREQ 

0.1075+J2.854 

1. 3604- j  13. 84 

1.401+331.93 

IkHz 

1.284  +j26.37 

10.77+jll5.9 

23.11+3279.9 

lOkHz 

6.385  +j249.2 

41. 84+j 1060.0 

93.83+32556.0 

lOOkHz 

23.14  +J2440.0 

113.3+310298.0 

200.4+324854.0 

IMHz 

BENCHMABK  CALCULATION 

INDUCTANCE  OF  COIL  IN  AIR 

0.471  (-6)  H 

(OUR  MODEL) 

0.473  (-6)  H 

(GROVER,  "INDUCTANCE 

CALCULATIONS,"  VAN  NOSTRAND,  1946) 

ALL  CALCULATIONS 

ASSUME  A  CLOSE-PACKED 

COIL  OF  20  GAUGE  HIBE, 

CONSISTING  OF 

15.3  TURNS. 
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A  more  common  presentation  of  the  data  in  Table  1  la  the  iiiq»edance** 
plane  plot  of  Figure  9.  Here  the  real  and  Imaginary  parts  of  the  driving- 
point  Impedance  are  displayed,  with  frequency  as  a  parameter.  The  three 
curves  correspond  to  the  coll  alone,  coll+cyllndrlcal  core,  and  coll-f’E- 
shaped  core,  each  In  the  presence  of  the  wrkplece.  The  interpretation 
of  the  diagram  Is  that  the  coupling  to  the  workpiece  becoe^s  tlghtar  in 
going  from  the  coll  alone  to  the  coll+E-shaped  core. 

The  reactive  part  of  the  driving-point  impedance  Is  due  partly  to 
the  self-inductance  of  the  coil+core  combination  (as  If  the  workpiece 
were  absent)  and  partly  to  mutual  Interaction  of  the  Induced  eddy-currents 
within  the  workpiece.  The  resistive  part  of  the  Impedance  Is  due  entirely 
to  the  effects  of  the  Induced  eddy-currents.  At  higher  frequencies  the 
eddy-currents  reside  nearer  the  surface  of  the  workpiece,  due  to  the  skin- 
effect,  therefore  Tnflk^ng  the  volume  of  Interaction  within  the  workpiece 
smaller.  This  means  that  the  volume  losses  within  the  workpiece  decrease, 
thereby  decreasing  the  effective  resistance  of  the  workpiece,  making  the 
driving-point  impedance  of  the  coil+core  combination  reactive. 

(c)  Induced  Fields  Within  the  Workpiece 

The  purpose  of  the  probe  coil  Is  to  Induce  eddy-currents  within  the 
workpiece.  It  Is  a  fairly  straight-forward  computation  to  determine  the 
distribution  of  Induced  eddy-currents,  once  the  induced  magnetization 
within  the  core  la  known.  The  fundamental  relation  Is 

•  dE^ir^z)  ,  (67) 

where  the  subscript  ”2"  denotes  region  2,  the  workpiece,  and  a  Is  the 
conductivity  of  the  workpiece.  Hence,  our  first  job  Is  to  determine  an 
Integral  expression  for  A2(r,z)  In  terms  of  the  current  and  magnetization 
In  region  1,  where  the  coll  and  core  are  located.  The  details  follow  the 
method  of  Section  II. 

We  start  by  replacing  (6)(a),(b)  by 
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V  .  V  .  -  ko  §12  -  0 


(r,z)  In  1 
(r',z')  in  2 


(68) (a) 


n  Tf  ,  2  -  _6(r-r')6(z-z')  - 

V  z  V  X  G22  -  G22  ^ 


(t,z)  In  2 
(r',z')  in  2 


Next,  ve  oultiply  (5)  by  the  appropriate  Green's  function  from  (68) 
and  then  use  the  vector  Green's  identity  over  regions  1  and  2»  as  before, 
to  get 

2iniQ  II  J(r,z)Gj^2(’^»*5  r',z')rdrd2  +  2iTyQ  ||  ^  *  M(r,z)  •  r',z')rdrdz 

coil  core 

■  -2iT|*[A2(r,z)  X  Vx  Gj^2(’f»*»  r' ,z')  -  G^, 2(^.2*  r' ,z')  x  Vx  Aj^(r,z)]  |  •  a^rdr  (69)1 


A2(r',z') --2Tr  |*(A2(r,z)  X  7  X  G22(r,z;  r',z')  - 
0 

-  G22(r,z;  r',z')  x  7  x  A2(r,z)]  |  *  a  rdr 

z“0 


which  replace  (8)  and  (9),  respectively. 

The  boundary  conditions  on  G^2  and  5^2  at  z  >  0  are: 


a  X  G-_  •  a  X  G_« 
z  12  z  22 


•—  a  X 7 X G-_  «  —  a  x7xG, _ 
Uj^  z  22  Mq  z  12 


Arguing  as  before  in  Section  11,  can  show  that  the  integral  in 
(69)  (b)  is  surface  integral  on  the  right-hand  side  of 

(69) (a).  Therefore  (69) (a)  may  be  written 
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A^Cr'.z')  -  aTfUj^  II  J(r,z)Gj^2^’^»*»  r’,z’)rdrdz 
coil 


+  liJV,  II  ^  *  M(r,z)  •  G^^(,T,zi  r',z')rdrdz  .  (70) 


core 


There  exists  a  reciprocity  relation  between  ®2l* 


UiGiiCr.z;  r*,z')  -  u„G,,(r,z;  r’.z’)  , 


O^ll' 


(71) 


and  with  this  (70)  becomes 


^^(r'.z')  -  2hUq  II  G2j^(r',z';  r,z) J(r,z)rdrdz 
coll 


+  2miQ  II  *  M(r,z)  •  Gj,j^(r’,z’;  r,z)rdrdz  .  (72) 

core 


This  is  the  desired  form  because  we  have  previously  computed  in 
(34(b). 

As  usual,  we  transfer  the  curl  operation  to  the  Green’s  function  in 
the  second  integral  of  (72)  by  using  (15),  and  then  we  substitute 

M  -  B(-i - r)! 

\  V 


^2  *2 


A2(r’,z')  -  y^n^I  J  rdr  J  dz  J  ^  e  ^  Jj^(r-t) Jj^(r 


ri  0 


ff  *^o  r 

+  11  rdrdzCl  -  ~)B^(r,z)  1  Te  ^  ° 


core 


ff  ^0  r  2 

JJ  rdrdzd  -  -jJ^)B^(r,z)  J  Te  ^  J^(r'^)JQ(r£)tV  .  (73) 


core 


f 
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When  (37)  is  substituted  into  (73),  and  the  Integrals  with  respect 
to  (r,z)  and  (r',z’)  are  evaluated,  we  get  the  final  expression 


r  Hr  Mr  1)  , 

A2(r',z')  -  j  - Z - Jj^(r'Z)dZ 


+  (1  -  e" 

H  J  J  J  J 


(74) 


where 


,(r) 


1  .  (-)  (+)  ct  z’ 

~J  T(£)F^(zj  ^ZJ  )  — - ^ ^ -  eV  JJ(r'^)d^  (75)(a) 


.(z) 


^  ) e“l*’  J^(r  'l)dl 


(b) 


Using  (74)  we  have  computed  the  induced  E-field  in  air,  and  the  in¬ 
duced  eddy-currents  in  the  workpiece,  as  a  function  of  radial  position 
at  two  depths,  z  ■  0  and  z  -  -2.7mm,  as  shown  in  Figure  10.  The  latter 
depth  is  the  skin-depth  in  aluminum  (a  -  3.5  x  10^)  at  IkHz. 

In  Figure  11  we  show  £  in  air  at  these  two  depths,  the  three  curves 
corresponding  to  the  three  core/coll  combinations,  coll  alone,  coll-f 
cylinder,  and  coil+E-shape.  The  reduction  in  field  intensity  at  -2.7mm 
is  not  due  to  skin  effect  (the  workpiece  is  absent) ,  but  is  simply  due 
to  distance;  this  is  a  manifestation  of  "lift-off."  It  is  clear  from 
part  (a)  of  the  figure  that  the  coll  alone  produces  the  weakest  field 
and  that  the  E-shaped  core  produces  the  strongest  peak.  This  is  to  be 
expected.  A  more  interesting  result,  however,  is  that  the  E-shapcd  field 
falls  off  very  quickly,  and  becomes  the  weakest  of  the  three  fields  beyond 
the  radial  limit  of  the  core.  This  Is  due  to  the  fact  that  the  net  mag¬ 
netic  flux  that  is  linked  by  a  circular  path  beyond  this  radius  Is  prac¬ 
tically  zero,  because  virtually  all  of  the  flux  leaving  the  outer  leg 
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returns  through  the  middle  leg.  Hence,  we  conclude  that  the  E-shaped  core 
Is  effective  for  producing  an  Intense,  concentrated  field.  The  results 
of  Figure  11  are  Independent  of  frequency  over  the  range  IkHz  to  IMHz. 

The  presence  of  the  workpiece  enhances  the  concentrating  effect  of 
the  coll/core  combinations,  as  Figure  12  shows.  Here,  the  Induced  eddy- 
currents  are  shown  at  IkHz.  The  fields  now  are  complex,  so  that  the  real 
and  imaginary  parts  are  both  shown,  the  Imaginary  part  dominating.  If  the 
curves  of  Figure  11  are  compared  with  the  curves  of  the  Imaginary  parts 
(marked  "I")  of  Figure  12,  It  will  be  seen  that  the  latter  curves  fall 
off  more  rapidly  with  r,  for  all  three  combinations  of  coll  and  core. 

This  enhanced  concentration  Is  due  to  Lenz'  law,  idilch  states  that  in¬ 
duced  currents  flow  In  such  a  direction  as  to  oppose  the  flux  that  pro¬ 
duces  them.  If  the  flux  Is  reduced  then  the  electric  field  and  eddy- 
currents  must  fall  off  more  rapidly  with  radial  distance.  The  reduction 
of  the  fields  at  -2.7mm,  as  shown  in  part  (b)  of  the  figure,  is  now 
principally  due  to  skin  effect. 

The  effects  of  increasing  frequency  are  shown  In  Figure  13,  which 
illustrates  the  real  and  imaginary  parts  of  the  induced  eddy-currcaits 
at  INHz.  At  this  frequency  these  parts  are  approximately  equal,  and 
the  fields  vanish,  due  to  skin-effect,  at  -2.7mm  (which,  it  will  be  re¬ 
called,  Is  the  skin-depth  at  IkHz).  The  rather  erratic  nature  of  the 
curves  that  correspond  to  the  coil+cyllnder  and  coil+E-shape  combinations 
is  probably  due  to  numerical  instabilities  in  computing  the  coefficient 
Integrals  and  of  (75) .  We  are  presently  Investigating  methods 

for  smoothing  these  data. 

In  all  of  the  cases  that  are  illustrated  In  Figures  11-13,  it  is 
clear  that  the  E-shaped  core  is  most  effective  in  coupling  energy  Into 
the  workpiece,  as  was  first  suggested  In  the  discussion  of  Figure  9. 

This  Is  another  indication  of  the  Internal  consistency  of  our  model. 
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VI.  COMMENTS,  EXTENSIONS  AND  CONCLUSIONS 


» 


Our  objective  In  this  report  was  to  apply  modern  methods  of  computa¬ 
tional  electromagnetics  to  develop  a  model  of  eddy-current  probe  colls 
with  ferrite  cores,  which  would  provide  a  systematic  and  rational  basis 
for  the  design  and  characterization  of  such  colls.  Ihe  model  that  we 
have  demonstrated  meets  this  goal;  In  addition,  it  Is  physically  and  math¬ 
ematically  consistent. 

Our  approach  was  to  ttse  a  volume  integral  equation  and  the  method  of 
moments.  Other  approaches  are  available,  namely  finite  elements,  and  a 
surface  Integral  equation  plus  the  method  of  moments.  The  volume  inte¬ 
gral  equation  approach  resulted  In  a  simpler  model,  that  required  only 
a  single  Green's  fxinctlon  (for  the  whole  of  regions  1  and  2),  rather  than 
a  Green's  function  for  the  Interior  of  the  ferrite  core  and  a  separate 
Green's  function  for  the  exterior,  as  would  have  been  required  in  the 
other  two  approaches.  This  resulted  In  simpler  programming  for  the  com¬ 
puter.  It  is  usually  true,  however,  that  different  problems  favor  dif¬ 
ferent  methods,  and  one  must  be  familiar  with  all  approaches  In  order  to 
select  the  one  that  is  most  efficient. 

The  model  that  we  have  derived  Is  useful  for  the  analysis  of  probes 
with  linear.  Isotropic  ferrite  cores.  It  can  be  extended  to  apply  to 
the  following  problem  areas,  all  of  which  are  iiiq>ortant  in  contemporary 
eddy-current  NDE  technology: 

*  Analyze  probe  coils  with  anisotropic  (single 
crystal)  cores 

*  Time-dependent  (pulsed)  eddy-current  problems 
(linear  core) 

*  Time-dependent  (pulsed)  eddy-current  problems 
(nonlinear  core) 

*  Optimization  of  probe  coll  design  to  satisfy 
simultaneous  constraints  on  eddy-current 
distribution  and  driving-point  impedance 

In  principle,  the  only  extension  that  Is  required  of  the  present 
model  In  order  to  treat  the  first  and  third  problems  Is  to  define  a 
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different  B  vs.  M  relationship  than  the  one  used  here.  For  the  single 
crystal  problem  this  relationship  would  remain  linear,  but  would  Involve 
a  tensor  permeability,  rather  than  the  scalar  p. 

The  nonlinear  problem  probably  can  be  approached  through  the  theory 
of  micromagnet Ics  [12],  which  would  require  the  simultaneous  solution  of 
a  volume  Integral  equation  of  the  type  described  In  this  report  and  a  non¬ 
linear  partial  differential  equation  (the  fundamental  equation  of  mlcro- 
magnetlcs) .  It  Is  possible  that  the  best  approach  to  the  numerical  com¬ 
putation  with  such  a  model  would  be  a  combination  of  finite  elements  and 
the  method  of  moments. 

All  of  these  problems  are  solvable;  they  only  require  a  sustained 
modeling  effort.  The  rewards  for  such  an  effort  could  be  immense  In  terms 
of  Improved  design  and  understanding  of  probe  colls  with  ferrite  cores. 
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Figure  1.  Ferrite  core  in  the  shape  of  a  general  body-of-revolution 
with  its  exciting  coil,  in  the  presence  of  a  stratified 
ha1f*space  workpiece. 
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Figure  2.  E-shaped  body-of-revolution  with  exciting 


44 


WORKPIECE 


f 

Figure  3.  Truncated  cylinder  with  exciting  coil. 
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REGION  2 


Figure  4.  Showing  the  two  regions  and  surfaces  that  are  used  In  the 

application  of  the  vector  Green's  Identity  to  the  derivation 
of  the  volume  Integral  equation. 
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Figure  5 


A  partitioning  of  the  cylindrical  core  into  a  regular 
grid,  for  applying  the  method  of  moments. 
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Figure  6.  Definition  of  1th  cell  boundaries,  for  computation  of 
Integrals. 


48 


Figure  8. 


Approximate  8-f1eld  within  the  E-shaped  core. 
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Workpiece  or  Air 
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Figure  10.  Showing  coll+core  and  region  of  workpiece,  or  air,  in 
which  eddy-currents,  or  E-field,  are  calculated. 
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Figure  11.  E-field  in  air  at  surface,  (a),  and  -2.7inm,  (b),  levels 
of  region  of  Figure  10. 
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Figure  12.  Induced  eddy-currents  In  workpiece,  at  1  kHz.  (a)  surface, 
(b)  -2.7tnni.  R»real  part,  I«1inag1nary  part.  Upper-left 
legend  Indicates  radial  extent  of  core-shapes  and  coll. 
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